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Abstract 

Wc present IntU package for Mathematica computer algebra system. The 
presented package performs a symbolic integration of polynomial functions 
over the unitary group with respect to unique normalized Haar measure. We 
describe a number of special cases which can be used to optimize the calcu- 
lation speed for some classes of integrals. We also provide some examples of 
usage of the presented package. 
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Nature of problem: Symbolic integration of polynomial functions over the unitary 
group with respect to unique normalized Haar measure. 

Solution method: A package of functions for Mathematica computer algebra sys- 
tem. 

Restrictions: Running time of the presented procedures grows rapidly with the 
dimensionality of the problem. 

Running time: For examples listed in the paper the running time is shorter than 
20 sec. 



1. Introduction 

The integration over unitary group is an important subject of studies 
in many areas of science, including mathematical physics, random matrix 
theory and algebraic combinatorics. In 2006 Collins and Sniady [IJ proved a 
formula for calculating monomial integrals with respect to the Haar measure 
on the Unitary group 

/ UjjUj,j,dU= I f/,,,,... E7,,,.;df/. (1) 

^U(d) JV{d) 

Integrals of the above type are known as moments of the U((i) and are well- 
known in mathematical physics literature for a long time. The problem of the 
integration of elements of unitary matrices was for the first time considered 
in the context of nuclear physics in [2]. The asymptotic behaviour of the 
integrals of the type ([T| was considered by Weingarten in [3]. 

In this paper we describe a Mathematica package IntU [1] for calculat- 
ing polynomial integrals over U(c?) with respect to the Haar measure. We 
describe a number of special cases which can be used to optimize the calcu- 
lation speed for some classes of integrals. We also provide some examples of 
usage of the presented package including the applications in the study of the 
geometry of the quantum states. 

This paper is organised as follows. In Section |2] we introduce notation 
present mathematical background concerning polynomial integrals over uni- 
tary group. In Section |3] we describe some special cases, in which the in- 
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tegration can be calculated more efficiently. In Section |4] we provide the 
description of the IntU package with the list of main functions. In Section [5] 
we show some examples of the usage. In Section [6] we provide a summary of 
the presented results and give conclusions. 

2. Mathematical background 

2.1. Basic concepts 

We denote by square matrices of size n. The compact group of dx d 
unitary matrices we denote as \J{d). We equip the above group with unique 
normalized Haar measure denoted by dU. Random elements distributed with 
measure dU form so called Circular Unitary Ensemble. 

Integer partition A of a positive integer n is a weakly decreasing sequence 
A = (Ai, A2, . . . , A/) of positive integers, such that Yl\=i^i — \M = ^- To 
denote that A is a partition of n we write A h n. The length of a partition 
is denoted by /(A). By A U yU we denote a partition of ni + n2 obtained by 
joining partitions A h ni and fi\- n2. 

Each permutation a E Sn can be uniquely decomposed into a sum of 
disjoint cycles where the lengths of the cycles sum up to n. Thus the vector 
of the lengths of the cycles, after reordering, forms a partition A h n. The 
partition A is called the cycle type of permutation a. 

2.2. Moments of the U((i) 

Let us consider a polynomial p. From the linearity of an integral we have 

/ p{U)dU= V c{I,J,I\f) [ UjjUrj.dU, (2) 

Jv{d) J\J(d) 

where /, J, /', J' are multi-indices and c are the coefficients of p. The value 
of such monomial integrals is given as [1] 

/ UuUrrdU = [ U,,j,...U,^,JJ^,^...m^,^dU= (3) 

Jv{d) Jv{d) 

= > 5i, i' ...5i i' Sj, ...5i i' Wq(Ta'^,d), 
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where Wg is the Weingarten function discussed below. In the case where the 
multi-indices differ in length, i.e. n ^ n', we have 

/ UijUrj.dU= [ Ui,,,...Ui„,„m>,>...U,^rdU = 0. (4) 

Jvid) J\J{d) " " 

The integrals of the above type are known as moments of the U((i). 

2. 3. Weingarten function 

The Weingarten function [T] is defined for a E Sn and positive integer c?, 

as 

>^«(''.'')=(;^E J^xV), (5) 



Ahn 
l{X)<d 



where the sum is taken over all integer partitions of n with length /(A) < d, 
SA,d(l) is the Schur polynomial sx evaluated at (1,1,...,!) and is an 

d 

irreducible character of the symmetric group Sn indexed by partition A. 

2.3.1. The dimension of irreducible representation of\]{d) 

The value of the Schur polynomial at the point (1,1,...,!), i.e. the 

d 

dimension of irreducible representation of U((i) corresponding to partition A, 
is equal to (see e.g. [HI Theorem 6.3]) 

^a,.(1) = .a(1,1,...,1)= n (6) 

d l<i<j<d 

2.3.2. Irreducible character of Sn 

The irreducible character of Sn indexed by partition A, X^i.^) depends 
on a conjugacy class of permutation a. Two permutations are in the same 
conjugacy class if and only if they have the same cycle type, thus it is common 
to write x^i.^) = X^{lA-i where fi is an integer partition corresponding to the 
cycle type of a. 

In the case of identity permutation the cycle type is given by a trivial 
partition, e = {l,l,...,l}, and the character value is equal to the dimension 

V 

n 

of the irreducible representation of Sn indexed by A. In this case it is given 
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by the celebrated hook length formula [5| Eq. 4.12] 



1 Li,j "-ij 



where |A| = Ai + A2 + ■ ■ ■ + ^i{x), and hfj is the hook length of the cell 
in a Ferrers diagram corresponding to partition A, see e.g. [6l p. 57]. 

In the case of a non-trivial partition the character of symmetric group 
X'^(cr) = can be evaluated with the use of Murnaghan-Nakayama rule 

(see e.g. [Zl Th. 4.10.2]), which describes a combinatorial way of calculating 
the character. [8] is used. 

From the above considerations one can notice that the Weingarten func- 
tion depends only on a cycle type of a permutation a and thus it is constant 
on a conjugacy class represented by a. Thus we may define the Weingarten 
function as 

;{A)<d 

where 11 is an integer partition, which is a cycle type of a. 



3. Special cases 

In this section we present some special cases of integrals, with respect 
to the Haar measure on the unitary group. In these cases the value of the 
integral can be calculated without the direct usage of formula (|3]), which 
requires processing of Ilf ^^'Ilj^i' permutations, where fcj (Ij) denotes the 
number of i (j) in multi-indices I (J) respectively. 

The presented special cases have been implemented in the package, to 
increase its efficiency. This goal has been achieved by minimizing the number 
of calculations of the Weingarten function. 
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3.1. First two moments ofl]{d) 

In the case of monomials of rank equal to 2 and 4, we have the following 
formulas P Prop. 4.2.3] 

UijUi'j'dU = -6ii'6jf, (9) 

U(d) " 

Ui^j^Ui^j^Ui'j'^Ui'^j^dU = — — (lUj 

u(d) a — i 



3232 



d (rf2 - 1) 

allowing us to calculate the values of polynomial integrals of degree less than 5 
without the direct usage of the formula ([3]). 

This optimization gives only a minor improvement of efficiency, as in these 
cases the direct calculation of Weingarten function is very fast. 

3.2. Elements from one row (column) 

The next optimization is based on the fact that the distribution of random 
vector consisting of squares of absolute values of a row (or a column) of 
unitary matrix distributed with the Haar measure, is uniform on a standard 
d-simplex A'=' [9, Ch. 4] 

{Ki|MM,,2|',...,|MMr}~u(A'^), (11) 

where U(A) denotes the normalized uniform measure (proportional to Lebesgue 
measure) on a set A C M'', and standard rf-simplex A*^ is defined as. A'' = 

Using Beta integral, one obtains that for a fixed row and a vector p 
with non-negative entries pj, we have the following 



V(d) 
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The above, as a special case, gives us: 



{d-l)\k\ 



U(d) 



u..rdU ^ (^3Tf^' (13) 



M'KM'dU = (14) 



U(d) 



which can be found in hterature P, HU] . 

This optimization allows for an enormous improvement in efficiency thanks 
to avoiding (^Pi)! IIp*' executions of Weingarten function needed in the case 
of the direct usage of formula ^ . 

3.3. Even powers of diagonal element absolute values 
In this subsection we consider the integrals of the type 

/ \u,j\''P\uk,i\'"'dU, (15) 
Jvid) 

where p, q are non-negative integers and i ^ k, j ^ I. In the case of p = g = 1 
this integral is known [9l Prop. 4.2.3] 

/ \uij\^\uk,i\^dU = (16) 
Jv{d) a - i 

For general non-negative integers p, q the following proposition is true. 

Proposition 1. Let p,q be non-negative integers, for i ^ k and j ^ I, we 

have 

[ Kjl'^M^'dU = VkAk^-W?(A U /x,ci), (17) 

where the above sum is taken over all integer partitions of p and q. Symbol 
ki, denotes a cardinality of conjugacy class for a permutation with cycle type 
given by partition u \- r 1^ Eg. 1.2] 

(18) 



l'^imi!2™2m2! ■ ■ ■ r^-m^!' 
where denotes the number of i in partition v. 
The above is a special case of more general fact. 
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Proposition 2. For any permutation vr G 5"^ and any non-negative integers 
Pi,P2, ■■ - Pd, we have 

„ d „ d 

/ ni%-o)i''^'^^= / ni«^-.^f'^'^f^ (19) 

^j=l / Aihpi A2I-P2 Adhpd \j=l / 

Proof. If we apply the formula from Eq. ^ to integral (19), then the non- 
vanishing permutations of indices are these which permute within the blocks 
of sizes {pi,p2, ■ ■ ■ ,Pd}, i-e. 

where aj G Sp. and it permutes indices in j*^ block of size pj. The same 
situation holds in the case of the second indices. Thus the permutation rcr~^ 
is also in this form. Each permutation of the above type will be present in 
the sum 11^'=! Pi' times. The cycle type of such permutations is given by 
a partition which is obtained by joining cycle types of small permutations. 
Since Weingarten function depends only on a cycle type of permutation, the 
size of conjugacy class is calculated for each partition, instead of evaluating 
multiple times Weingarten function. □ 



The formula (19) is far less computationally demanding than the direct 
usage of ([3]). One can notice that Proposition |2] allows us to avoid at least 
n^=iPi' executions of Weingarten function comparing to the direct usage of 
(|3|. However, the exact time-efficiency depends also on the cardinality of 
conjugacy classes for partitions of pi, ■ ■ ■ ,Pd- 

3.4- Cycle permutations 

In this section we consider another special type of integral for positive 
integer k and a permutation a of {1, 2, ... , m} being a cycle 

/ In^M^MO) dU. (20) 
Jv(d) \X\ I 



We have the following proposition. 
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Proposition 3. Let m < d, k is a positive integer and a G Sm be a cycle, 
i.e. the cycle type is given by partition {m}. Then 



TT I dU = {k\Y"'-'y2kxWg{m\,d). (21) 

i=l J Ahfc 

Proof. One can notice tliat cycle lengtlis of permutations in this setting must 
be divisible by m. The number of permutations with cycle type given by a 
particular partition can be easily obtained by the usual counting argument. 

□ 

Using the above formula ( [2l| ) one obtains an enormous improvement in 
efficiency by avoiding more than (/c!)^™"^ executions of Weingarten function, 
comparing to the case of the direct usage of ([3]). 

4. Package description 

Below we describe the functions implemented in IntU package. The func- 
tions are grouped in three categories: functions implementing the main func- 
tionality, functions related to the calculation of Weingarten function and 
helper functions. 

4.1. Main functionality 

The main functionalty of IntU package is provided be the IntegrateUnitaryHaar 
and IntegrateUnitaryHaarlndices functions. The first one operates di- 
rectly on polynomial expressions, while the second one accepts four-tuple of 
indices. The examples of the usage are given in Section [5j 

• IntegrateUnitaryHaar [integrand, {var, dim}] - gives the definite 
integral on unitary group with respect to the Haar measure, accepting 
the following arguments 

— integrand - the polynomial type expression of variable var with 
indices placed as subscripts, can contain any other symbolic ex- 
pression of other variables, 

— var - the symbol of variable for integration, 

— dim - the dimension of a unitary group, must be a positive integer. 



This function is presented in the examples described in Sections 5.1 
and lSTSl 
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IntegrateUnitaryHaar [f , {u,dl} , {v,d2} , . . .] - gives multiple in- 
tegral 

^ dU [ dV...f. (22) 

U(dl) iu(d2) 



This function is presented in the example described in Section 5.4 



• IntegrateUnitaryHaarlndices [{II , Jl , 12, J2} ,dim] - calculates the 
integral in Eq. ([s]) for given multi-indices II , Jl , 12, J2 and the dimen- 
sion dim of the unitary group. This function is presented in the example 
described in Section 15.21 

4-2. Weingarten function 

The main functions implemented in the package, IntegrateUnitaryHaar 
and IntegrateUnitaryHaarlndices, utilize the following functions to find 
the value of the integral. 

• Weingarten [type , dim] - returns the value of the Weingarten function 
given in Eq. (|8|) and accepts the following arguments 



type - an integer partition which corresponds to cycle type of 



permutation (see Section 2.3) 



— dim - the dimension of a unitary group, must be a positive integer. 
CharacterSymmetricGroup [part , type] - gives the character of the 



symmetric group ^^^''^(type) (see Section 2.3.2). 

Parameter type is optional. The default value is set to a trivial par- 
tition and in this case the function returns the dimension of the irre- 
ducible representation of symmetric group indexed by part, given by 
Eq. ([T]). If type is specified the value of the character is calculated by 
Murnaghan-Nakayama rule using MNInner algorithm provided in [8]. 

SchurPolynomialAtl [part , dim] - returns the value of the Schur poly- 
nomial Spart at point (1,1,..., 1), i.e. the dimension of irreducible 

dim 

representation of U(dim) corresponding to part, see Eq. ([6]). 
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4 ■ 3. Helper functions 

• PermutationTypePartition [perm] - gives the partition which repre- 



sents the cycle type of the permutation perm (see Section 2.1). 



MultinomialBeta [p] - for a given ci-dimensional vector of non-negative 
numbers Pi,P2, ■ ■ ■ ,Pd returns the value of multinomial Beta function 
defined as 

r(Ei=iPi) 



This function is used in the optimization described in Section 3.2 



• ConjugatePartition [part] - gives a conjugate of a partition part 
(see ^). This function is used for calculating hook- length formula 
given by Eq. ([T]). 

• CardinalityConjugacyClassPartition [part] - gives a cardinality 
of a conjugacy class for the permutation with cycle type given by par- 
tition part (see [71, Eq. 1.2]). This function is used in the optimization 
described in Section 13.31 

• BinaryPartition [part] - gives a binary representation of a parti- 
tion part. This function is needed for the implementation of MNInner 
algorithm. 

5. Usage examples 

In order to present the main features of the described package we provide 
a series of examples. 

5.1. Elementary integrals 

Let us assume that d = 3. We want to calculate the following integrals 

ui,i\^dU, (24) 

V{d) 

M'M^dU, (25) 

U(d) 

f 

U(d) 



/ Ul,lU2,2Ul,2U2,ldU. (26) 
Jv(d) 
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Let us start by initializing the package 

ln[1]:= « IntP' 

Package IntU version 0.2.0 {last modification: 19/09/2011). 

ow we calculate the integrals. 

In[21:= d = 3; 

in[3i:= IntegrateUnitaryHaar [Abs[ui ^ 2 , {u , d}] 
1 

Clut[3]= - 

3 

in[4]:= IntegrateUnitaryHaar ^Abs^Uj^ 2 "2,2] ^ r ("^ 
1 

Clut[4]= - 

8 

in[5i:= IntegrateOnitaryHaar [ Uj I ^ Conjugate [u^ 2 "2,1]/ 
1 

Qjt[5]= 

24 



5.2. Operations on indices 

Let us take the following set of multi-indices 

J = {LLh2,2}, J = {2,2, 1, 1, 1} (27) 
/' = {1,1,1,2,2}, J' = {2, 1,1, 2,1} (28) 

and set d = 6. The above is equivalent to expression 

'"l,2Ml,2Ml,lM2,lM2,lMl,2Ml,lMl,lM2,2M2,l (29) 

with symbolic variable which we aim to integrate over U(d). After simpli- 
fication the expression is equal to 

|Ml,inMl,2n'"2,lPwi,2W2,lMl,l«2,2- (30) 

After initializing the package and defining appropriate indices 

ln[2]:= Jl= {I, 1, 1, 2, 2); Jl= {2, 2, 1, 1, 1); 
12= {1, 1, 1, 2, 2}; J2= {2, 1, 1, 2, 1}; 
d= 6; 

we calculate the integral using IntegrateUnitaryHaarlndices as 
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in[5]:= IntegrateUnitaryHaarIndices[{Il, Jl , 12, J2), d] 
I 

Qjt[5]= 

16200 

which is equivalent to 

in[6]:= IntegrateUnitaryHaar [ 

Abs "1,2 ''2,1]^^ ^1,2 ^2,1 Conjugate ["2,1 '^2, 2]^ {"/ <^}] 

I 

Qjt[6]= 

16200 



5.3. Matrix expressions 

IntU package allows us to integrate matrix expressions, for example let 
us take d = 2 and integrate 



(31) 



U(d) 



ln[2] 



We define symbolic matrices U G \]{d) and U2 = U ® U G V{d^) as 

d = 2; 

C7 = Array [ u„j „2 4, {d, d}]; 
[72 = KroneciterProduct [7, P]; 

and construct the integrand as 

in[5]:= integrand = KroneckerProduct[lI2 , Conjugate [ C72 ]] ; 

By using IntegrateUnitaryHaar function 

in[6]:= IntegrateUnitaryHaar [integrand, {u, d}] 



we learn that the integral in Eq. (31) is equal to 
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(32) 
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5.4- Mean value of local unitary orbit 

In this example we calculate the mean value of a local unitary orbit of 
a given matrix X G A^d2 

E[{U ®V)X{U ®V)^]. (33) 

We assume that U and V are stochastically independent random unitary 
matrices of size d distributed with the Haar measure. In this case we have 

E[{U^V)X{U®Vy]= [ [ {U ^V)X{U ^VydUdV. (34) 

In this example we take d = 3. 

We define symbolic matrices X of size d^ and U,V & l]{d) as 

ln[2]:= d = 3; 

X = Array [x„j„2 i, {d^2, d -2}]; 
U = Array [u„i „2 {d, d}]; 
V = Array [v3j „2 i, {d, d}]/ 

Using IntegrateUnitaryHaar function with two variable specifications we 
calculate the double integral 

in[6]:= int = IntegrateUnitaryHaar [ 

KroneckerProduct [U , V].X . KroneckerProduct [U , V ] * 
, {u, d}, {V, d}] 



and find out that the expectation value in Eq. (33) is equal to 



E[(f/ ® V)X{U ® VY] = ^tiX Idxd. (35) 

Similarly one can calculate the covariance tensor of the local unitary orbit 
of a given matrix X G Md2. U z = {U ® V)X{U ® VY then the covariance 
tensor is given by jTl] 

^{z^jZkl}i,kl] = HiU ® V)X{U ® Vy ® {U^V)X{U^Vy]. (36) 

Using IntegrateUnitaryHaar one can check that for the fixed dimension 
the integral agrees with the calculations presented in [TT] . 



5. 5. Moments of maximally entangled numerical shadow 

If [/ is a random dxd unitary matrix distributed according to the Haar 
measure then the pure state obtained by vectorization |,^) = -^\ec{U) is 
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maximally entangled on "H^ ®1-Lb = C'^ C'^. Moreover, state |^) has a 
distribution invariant to multiplication by local unitary matrices. We denote 
the corresponding probability measure on pure states of size d x d a.s fi. 
The numerical shadow of operator X G A^d2 with respect to this measure 
(maximally entangled numerical shadow) is defined as 

Px{z)= [ 5(z-(^|X|^))d/i(^). (37) 

The corresponding probability measure is denoted by df^x. For definition 
and basic facts concerning numerical shadows see [HI [T21 [13] • 

The first two moments of dfix are calculated in [TT], and are given by 

j^zdiij,iz) = ^trX , (38) 



and 



(39) 



l^zzdf^'xi^) = (|trXp+||X||y 

- ^3(^2^1) (l|t^^WllHS+l|trB(X)||y, 

where tr^ and tr^ denote partial traces over a specified sub-system and || • ||hs 
is a Hilbert-Schmidt norm given by ||X||hs = VtrXXt. 



In order to calculate (39) and (39) define symbolic matrices 

ln[2]:= d = 3; 

X = Array[x„,^„2 {d -^2, d ^2}]; 

17 = Array [ Uaj „2 i, {d, d}]; 

S = 1/ Sqrt I d ] Flatten [U]; 

z = S ■ X . Conjugate [S] ; 

zz = Simplify [ z Conjugate [z]]; 

Now we calculate the first moment 

in[8]:= IntegrateUnitaryHaar [ z , {u , d}] 

and the second moment 

in[9]:= IntegrateUnitaryHaar [zz , {u , d}] 

After some algebraic manipulations one can see that the above agrees with 



Eqs. (38) and (39). 
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6. Summary 



We described IntU package for Mathematica computing system for cal- 
culating polynomial integrals over U{d) with respect to Haar measure. We 
described some number of special cases which can be used to optimize the 
calculation speed for some classes of integrals. We also provide some exam- 
ples of the usage of the presented package including the applications in the 
geometry of the quantum states. 

Calculation time of the package strongly depends on a degree of the inte- 
grand. For polynomials of small degree the package is able to calculate the 
value of integral using the direct formula For polynomials of large degree 
the calculation time grows rapidly and the calculation is possible only if one 
of the special cases (optimizations) is used. 

Nevertheless, the presented package can be very useful in the investi- 
gations involving Circular Unitary Ensemble and the geometry of quantum 
states. 
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